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In this work we aim to highlight a close analogy between cooperative behaviors in chemical kinetics and 
cybernetics; this is realized by using a common language for their description, that is mean-field statistical 
mechanics. First, we perform a one-to-one mapping between paradigmatic behaviors in chemical kinetics 
(i.e., non-cooperative, cooperative, ultra-sensitive, anti-cooperative) and in mean-field statistical mechanics 
(i.e., paramagnetic, high and low temperature ferromagnetic, anti-ferromagnetic). Interestingly, the 
statistical mechanics approach allows a unified, broad theory for all scenarios and, in particular, 
Michaelis-Menten, Hill and Adair equations are consistently recovered. This framework is then tested 
against experimental biological data with an overall excellent agreement. One step forward, we consistently 
read the whole mapping from a cybernetic perspective, highlighting deep structural analogies between the 
above-mentioned kinetics and fundamental bricks in electronics (i.e. operational amplifiers, flashes, 
flip-flops), so to build a clear bridge linking biochemical kinetics and cybernetics. 



c 



ooperativity is one of the most important properties of molecular interactions in biological systems and it 
is often invoked to account for collective features in binding phenomena. 



In order to investigate and predict the effects of cooperativity, chemical kinetics proved to be a fundamental 
tool and, also due to its broadness over several fields of biosciences, a number of cooperativity quantifiers (e.g. 
Hills number 1 , Koshland cooperativity test 2 , global dissociation quotient 3 , weak and strong fine tunings 4 , etc.), 
apparently independent or distinct, have been introduced. However, a clear, unified, theoretical scheme where all 
cooperative behaviors can be framed would be of great importance, especially in biotechnology research 5,6 and for 
scientists dealing with interdisciplinary applications 7,8 . To this task, statistical mechanics offers a valuable 
approach as, from its basic principles, it aims to figure out collective phenomena, possibly overlooking the details 
of the interactions to focus on the very key features. Indeed, a statistical mechanics description of reaction kinetics 
has already been paved through theoretical models based on linear Ising chains 9 , spin lattices with nearest 
neighbors interactions 10 , transfer matrix theory' 10 and structural probabilistic approaches 11 . 

In this work we expand such statistical mechanics picture toward a mean-field perspective 12 by assuming that 
the interactions among the system constituents are not limited by any topological or spatial constraint, but are 
implicitly taken to be long-ranged, as in a system that remains spatially homogeneous. This approach is naturally 
consistent with the rate-equation picture, typical of chemical kinetics investigations and whose validity is restricted 
to the case of vanishing correlations 13,14 and requires a sufficiently high spatial dimension or the presence of an 
effective mixing mechanism (hence, ultimately, long-range interactions). In general, in the mean-field limit, 
fluctuations naturally decouple from the volume-averaged quantities and can be treated as negligible noise. 

By adopting a mean-field approach, we abandon a direct spatial representation of binding structures and we 
introduce a renormalization of the effective couplings. The reward lies in a resulting unique model exhibiting a 
rich phenomenology (e.g. phase transitions), which low-dimensional models typically lack, yet being still feasible 
for an exact solution. In particular, we obtain an analytical expression for the saturation function which is 
successfully compared with recent experimental findings, taken from different (biological) contexts to check 
robustness. Furthermore, from this theory basic chemical kinetics equations (e.g. Michaelis-Menten, Hill and 
Adair equations) are recovered as special cases. 

Further, there is a deep theoretical motivation underlying the development of a mean-field statistical 
mechanics approach to chemical kinetics: it can be used to code collective behavior of biosystems into a cyber- 
netical framework. In fact, cybernetics, meant as the science dedicated to the understanding of self-organization 
and emergent communication among the constituents of a system, can be naturally described via (mean-field) 
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statistical mechanics 15-17 . Thus, the latter provides a shared form- 
alism which allows to automatically translate chemical kinetics into 
cybernetics and vice versa. In this perspective, beyond theoretical 
interest, at least two concrete benefits may stem from our investiga- 
tion: first, in the field of biotechnologies, logical gates have already 
been obtained through biological hardware (see e.g. 5,6 ) and for their 
proper functioning signal amplification turns out to be crucial. In 
this paper, cooperativity in reaction kinetics is mapped into amp- 
lification in electronics, hence offering a solid scaffold for biological 
amplification theory. 

Then, as statistical mechanics has been successfully applied in the 
world of computing (for instance in neural networks 18 , machine 
learning 19 or complex satisfiability 20 ), its presence in the theory of 
biological processors could be of relevant interest. In particular, we 
discuss how to map ultra-sensitive kinetics to logical switches and 
how to read anti-cooperative kinetics as the basic ingredient for 
memory storage in biological flip-flops, whose interest resides in 
several biological machineries as gene regulatory networks 21 , ribos- 
witches 22 , synaptic switches 23 , autopoietic systems 24 and more 25-27 . 

To summarize, a rigorous, promising link between cybernetics and 
collective biological systems can be established via statistical 
mechanics and this point will be sketched and corroborated by 
means of several examples throughout this paper, which is structured 
as follows: 

First, we review the main concepts, facts and methods from both 
chemical kinetics and statistical mechanics perspectives. Then, we 
develop a proper theoretical framework able to bridge statistical 
mechanics and chemical kinetics; the former can also serve as a 
proper tool for describing and investigating cybernetics, thus, as a 
syllogism, chemical kinetics and cybernetics become also related. The 
agreement of our framework with real data, carefully extrapolated 
from recent biological researches, covering the various standard 
behaviors in chemical kinetics, is also successfully checked. Finally, 
results and outlooks are discussed. 



constant K w =K^\/ and dissociation constant K [i '=K^[ 
Kv{. Focusing on the steady state we get, iteratively, 



AD 



[pj-m 



Unfortunately, measuring [Pj] is not an easy task and one usually 
introduces, as a convenient experimental observable, the average 
number S of substrates bound to the protein as 



S = 



[P l ]+2[P 2 ]+ ...+N[P, 



Nj 



[Po] + [*]+. ■ 
K^[S]+2K^[S] 



Pn] 

2 , i it rWMUrlN 



-»#)[S]' 



(2) 



l+K( 1 )[S]+K( l )K( 2 )[S] 2 + ... +K( 1 )K( N )[SY 



which is the well-known Adair equation 3 , whose normalized express- 
ion defines the saturation function Y = S/N. In a non-cooperative 
system, one expects independent and identical binding sites, whose 
steady states can be written as (explicitly only for; = 1 and j = 2 for 
simplicity) 

0=-NK + [P 0 ][S]-K_[P 1 ], (/=!), (3) 



0 = 



-(N-l)K+[P 1 ][S]+2K_[P 2 ], (/=2), 



(4) 



where K + and K_ are the rates for binding and unbinding on any 
arbitrary site. Being K = K + /K_ the intrinsic association constant, 
we get 

[Pi] _ 
N[P 0 ][S] N ' 



K- 



0 = i), 



(5) 



Results 

Collective behaviors in chemical kinetics. Many polymers and 
proteins exhibit cooperativity, meaning that their ligands bind in a 
non-independent way: if, upon a ligand binding, the probability of 
further binding (by other ligands) is enhanced, like in the 
paradigmatic case of hemoglobin 9 , the cooperativity is said to be 
positive, vice versa there is negative cooperativity when the 
binding of more ligands is inhibited 28 , as for instance in some 
insulin receptors 29,30 and most G-protein coupled receptors 23,31 . 
Several mechanisms can be responsible for this effect: for example, 
if two neighbor docking sites on a polymer can bind charged ions, the 
electrostatic attraction/repulsion may be the cause of a positive/ 
negative cooperativity. However, the most common case is that the 
binding of a ligand somehow modifies the structure of the hosting 
molecule, influencing the binding over the other sites and this is the 
so-called allosteric mechanism 32 . 

Let us now formalize such behavior by considering a hosting 
molecule P that can bind N identical molecules S on its structure; 
calling Pj the complex of a molecule P with j e [0,N] molecules 
attached, the reactions leading to the chemical equilibrium are the 
following 

S+Pj-^Pj, 

hence the time evolution of the concentration of the unbounded 
protein P 0 is ruled by 

d[P 0 ] 



dt 



-K^[P 0 ][S]+K^[ Pl ], 



(1) 



where K+ v K^l\ are, respectively, the forward and backward rate 
constants for the state; = 1, and their ratios define the association 



K= 



2\P 2 



2K® 



(N-1)[P,][S] {N-lY 



0=2), 



(6) 



and, in general, = (N — j + l)K/j. Plugging this expression into 
the Adair equation (2) we get 



S= NK[S \ ^ Y KlS 



1+K\S] 



1+K\SV 



(7) 



which is the well-known Michaelis-Menten equation 3 . 

If interaction among binding sites is expected, the kinetics 
becomes far less trivial. Let us first sketch the limit case where inter- 
mediates steps can be neglected, that is 

[P 0 ]+N[S]^±[P N ], 



then 



S = 



N\P K 



NK\SY 



[Po] + [Pn] l + fSl J 



Y = 



K[Sf 



N l- 



(8) 



(9) 



More generally, one can allow for a degree of sequentiality and write 

K\S] n " 



Y = 



l + \S}' 



(10) 



which is the well-known Hill equation 3 , where n H , referred to as Hill 
coefficient, represents the effective number of substrates which are 
interacting, such that for n H = 1 the system is said to be non-cooperative 
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and the Michaelis-Menten law is recovered; for n H > 1 it is coopera- 
tive; for )if[^>l it is ultra sensitive; for n H < 1 it is anti cooperative. 

From a practical point ofview, from experimental data for T([S]), 
one measures n H as the slope of log(Y/(l — Y)) versus [S]. 

Mean-field statistical mechanics. One of the best known statistical 
mechanics model is the mean-field Ising model, namely the Curie- 
Weiss model 33 . It describes the macroscopic behavior of a magnetic 
system microscopically represented by N binary spins, labeled by i = 
1, 2, N, and whose state is denoted by a, ■= ±1. In the presence of 
an external field h and being J the N X N symmetric matrix encoding 
for pairwise interactions among spins, the (extensive, macroscopic) 
internal energy associated to a the configuration {a} = {o^, ct 2 > •■■» 
a N } is defined as 

£({a}|J, h)= - ^$><^- (11) 



It is easy to see that the spin configurations leading to a lower energy 
are those where spins are aligned with the pertaining field, i.e. Cj/z; > 
0, and pairs (i, j) associated to positive (negative) coupling are 
parallel (anti-parallel), i.e. a t aj = 1 (ff,ffy = —1). Notice that, in eq. 
11, we implicitly assumed that any arbitrary spin possibly interacts 
with any other. This is a signature of the mean field approach which, 
basically, means that interactions among spins are long-range and/or 
that the time-scale of reactions is longer than the typical time for 
particles to diffuse, in such a way that each spin/particle actually sees 
any other. We stress that the mean field approximation also implies 
that the probability distribution P({ff}) for the whole configuration is 
factorized into the product of the distribution for each single 

-r-rN 

constituents, namely P({a}) =ll i=1 P(<T,-), analogously to classical 
chemical kinetic prescriptions 10 . 

For the analytic treatment of the system it is convenient to adopt a 
mesoscopic description where the phase space, made of all the 2 N 
possible distinct spin configurations, undergoes a coarse-graining 
and is divided into a collection T of sets, each representing a meso- 
scopic state of a given energy E^, (we dropped the dependence on the 
parameters J, h to lighten the notation). In this way, all the micro- 
scopic states belonging to the set k e T share the same value of 
energy E^, calculated according to (11). In order to describe the 
macroscopic behavior of the system through its microscopical 
degrees of freedom, we introduce a statistical ensemble 
p = {Pk}keT' mean t as the probability distribution over the sets in 
T; consequently, p { s 0 and pk = 1 must be fulfilled. 
Accordingly, the internal energy and the entropy read as 

E{p)=J2 Ek Pk> S{p)=-K B J2pk^ogp k , (12) 



fe'7 



keT 



where K B is the Boltzmann constant, hereafter set equal to 1. Being fi 
> 0 the absolute inverse temperature of the system, we define the 
free-energy 

F(p,(S)=E(p)-(r 1 S(p). (13) 

Notice that the minimum of F ensures, contemporary, the minimum 
for E and the maximum for S, hence it provides a definition for the 
thermodynamic equilibrium. As a consequence, from eq. 13 we cal- 
culate the derivative with respect to the probability distribution and 
require dFldpt = 0; the solution, referred to as p, reads as 

r-fiBk 



Pk = 



(14) 



and it is called the Maxwell-Bolzmann distribution. The normaliza- 
tion condition implies Z(f}) = e~P Et and this quantity is called 
"partition function". We therefore have 



p 



■■p{E,-p,p)=-r l ^z{p). 



(15) 



In general, given a function /({ff}), its thermal average is (/) = 

£ W /(MK WM) / Z ^- 

As this system is expected to display two different behaviors, an 
ordered one (at low temperature) and a disordered one (at high 

temperature), we introduce the magnetization m = J ^ Gj/N, 
which provides a primary description for the macroscopic behavior 
of the system. In particular, it works as the "order parameter" and it 
characterizes the onset of order at the phase transition between the 
two possible regimes. More precisely, as the parameters f}, J, h are 
tuned (here for simplicity Jy > 0, Vi_ j), the system can be either 
disordered (i.e. paramagnetic), where spins are randomly oriented 
and (m) = 0, or ordered (i.e. ferromagnetic), where spins are con- 
sistently aligned and (m) ^ 0. The phase transition, separating 
regions where one state prevails against the other, is a consequence 
of the collective microscopic interactions. 

In a uniform system where h = J, V, j, /„ = 0, and hj = h, Vi, all 
spins display the same expected value, i.e. (ct;) = (a), Vi, which also 
corresponds to the average magnetization (m). Remarkably, in this 
case the free energy of the system can be expressed through (m) by a 
straightforward calculation 18 that yields 



F{fi,H) = -jr l ln2- 



-/T 1 lncosh[/i(/(m) + /z)]+^(m 2 }, (16) 



whose extremization w.r.t. to (m) ensures again that thermodynamic 
principles hold and it reads off as 



(m(/y,fc))=tanh0(7(ffl)+fc)] 



(17) 



which is the celebrated Curie- Weiss self-consistency. By simply solv- 
ing eq. 17 (e.g. graphically or numerically) the macroscopic behavior 
can be inferred. Before proceeding, we fix p = 1, without loss of 
generality as it can be reabsorbed trivially by h —* fih = h and / — » 
/»/ = /. 

In the non-interacting case (/ = 0), eq. 17 gets m(J = 0, h) = 
tanh(?z), which reminds to an input-output relation for the system. 
When interactions are present (/ > 0), one can see that the solution of 
eq. 1 7 crucially depends on /. Of course, (m(J = 0, h)) < (m(J > 0, h)), 
due to cooperation among spins, and, more remarkably, there exists a 
critical value J c such that when / s J c the typical sigmoidal response 
encoded by (m(J, h)) possibly becomes a step function (a true dis- 
continuity is realized only in the thermodynamic limit N— * °°, while 
at finite N the curve gets severely steep but still continuous). Hence, 
to summarize, the Curie Weiss model exhibits two phases: A small- 
coupling phase where the system behaves paramagnetically and a 
strong- coupling phase where it behaves ferromagnetically. The fer- 
romagnetic states are two, characterized by positive and negative 
magnetization, according to the sign of the external field. 

Statistical mechanics and chemical kinetics. In this section we 
develop the first part of our formal bridge and show how the 
Curie- Weiss model can be looked from a biochemical perspective. 
We will start from the simplest case of independent sites and later, 
when dealing with interacting sites, we will properly generalize this 
model in order to consistently include both positive and negative 
cooperativity. 

The simplest framework: non interacting sites. Let us consider an 
ensemble of elements (e.g. identical macromolecules, homo-allos- 
teric enzymes, a catalyst surface), whose interacting sites are overall 
N and labelled as i = 1, 2, N. Each site can bind one smaller 
molecule (e.g. of a substrate) and we call a the concentration of the 
free molecules ([S] in standard chemical kinetics language). We 
associate to each site an Ising spin such that when the i"' site is 
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occupied a, ■ = + 1, while when it is empty oy = — 1. A configuration 
of the elements is then specified by the set {a}. 

First, we focus on non-collective systems, where no interaction 
between binding sites is present, while we model the interaction 
between the substrate and the binding site by an "external field" h 
meant as a proper measure for the concentration of free-ligand mole- 
cules, hence h = h{a). One can consider a microscopic interaction 
energy given by 



E(a,h) = -h^Oi. 



(18) 



Note that writing h 0"; instead of ^^.ftj<7; implicitly assumes 
homogeneity and thermalized reactions: each site displays the same 
coupling with the substrate and they all had the time to interact with 
the substrate. 

We can think at h as the chemical potential for the binding of 
substrate molecules on sites: When it is positive, molecules tend to 
bind to diminish energy, while when it is negative, bound molecules 
tend to leave occupied sites. The chemical potential can be expressed 
as the logarithm of the concentration of binding molecules and one 
can assume that the concentration is proportional to the ratio of the 
probability of having a site occupied with respect to that of having it 
empty. In this simple case, and in all mean-field approaches 33 , the 
probability of each configuration is the product of the single inde- 
pendent probabilities of each site to be occupied and, applying the 
Maxwell-Boltzmann distribution P{a t = ±1) = e ±h , one finds 



P(a,- 



, + h 



and we can pose 



■P(ff,= -1) 



h — - log i. 
2 s 



(19) 



(20) 



The mean occupation number (close to the magnetization in statist- 
ical physics) reads off as 



«(W)=~X>+".-) 



Therefore, using ergodicity to shift S into (S) (see eq. 7), the satura- 
tion function can be written as 



Y(«) 



<s(W)> 



N 



= ^£^(M)«(M) = ^(i + W>0>)- (2i) 



amount of parameters, namely /, a only.], where a ligand belonging 
to one group interacts in a mean field way with ligands in the other 
group (cooperatively or competitively depending on the sign of the 
coupling, see below), and they both interact with the substrate. As a 
result, given the parameters J and h, the energy associated to the 
configuration {a} turns out be 



] N A N B N A + N„ 



(23) 



Some remarks are in order now. First, we stress that in eq. 23 the 
sums run over all the binding sites. As we will deal with the ther- 
modynamic limit (N— > °°), this does not imply that we model macro- 
molecules of infinite length, which is somehow unrealistic. Rather, 
we consider N as the total number of binding sites, localized even on 
different macromolecules, and the underlying mean-field assump- 
tion implies that binding sites belonging to the same group are all 
equivalent, despite some may correspond to the bulk and others to 
the boundaries of the pertaining molecule; such differences can be 
reabsorbed in an effective renormalization of the couplings. In this 
way the system, as a whole, can exhibit (anti) cooperative effects, as 
for instance shown experimentally in 34 . 

Moreover, for the sake of clearness, in the following we will assume 
that couplings between sites belonging to different groups are all the 
same and equal to /and, similarly, h t = h, for any i. This homogeneity 
assumption allows to focus on the simplest cooperative effects and 
can be straightforwardly relaxed. 

We also notice that this two-groups model can mimic both coop- 
erative and non-cooperative systems but, while for the former case 
bipartition is somehow redundant as qualitatively analogous results 
are obtained by adopting a fully- connected structure, for the latter 
case the underlying competitive interactions intrinsically require a 
bipartite structure. 

Now, the two groups are assumed as equally populated, i.e., N A = 
N B = N/2, such that their relative densities are p A = N A /N = p B = 
N B /N = 1/2. The order parameter can be trivially extended as 



m A -- 



(24) 



and, according to statistical mechanics prescriptions, we minimize 
the free energy coupled to the cost function (23) and we get, in the 
thermodynamic limit, the following self- consistencies 



(m A ) =tanh[/p„(m B ) + h], 



(25) 



By using eq. 17 we get 



Y(a) = -[l+tanh(fc)] 



(22) 



which, substituting 2h = log a, recovers the Michaelis-Menten beha- 
vior, consistently with the assumption of no interaction among bind- 
ing sites (/ = 0) in eq. 18. 

A refined framework: two-sites interactions. We now focus on pair- 
wise interactions and, seeking for a general scheme, we replace the 
fully- connected network of the original Curie- Weiss model by a 
complete bipartite graph: sites are divided in two groups, referred 
to as A and B, whose sizes are N A and N B (N = N A + N B ), respect- 
ively. Each site in A (B) is linked to all sites in B (A), but no link within 
the same group is present. With this structure we mirror dimeric 
interactions [Note: Notice that for the sake of clearness, we intro- 
duced the simplest bipartite structure, which naturally maps dimeric 
interactions, but one can straightforwardly generalize to the case of 
an ra-mer by an «-partite system and of course values of p A p B can 
be considered too. We did not perform these extensions because we 
wanted to recover the broader phenomenology with the smaller 



(m B ) = tanh[/p A (m A ) + h] 



(26) 



Through eqs. 25 and 26, the number of occupied sites can be com- 
puted as 

»A(W)=X)x(l + ff i )=N A : 



l + m A 



i=l 
N, 



n B 



(W)= Eki 



2 

l + m B 



(27) 



from which we get the overall binding isotherm 

(n A (a)) + (n B [tx)} 



Y(a)=- 



N 



(28) 



From eqs. 25-28 one can see that Y (a) fulfills the following free- 
energy minimum condition 



Y{a.-J) = -{ 1 + tanh 



/(2Y-l) + -log« 



(29) 
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Figure 1 | Theoretical predictions of typical binding isotherms obtained from the statistical mechanics approach (see eq. 28) are shown versus the 
substrate concentration a (main plot) and versus the logarithm of the concentration (inset). Different colors refer to different systems (hence different 
coupling strengths /), as explained by the legend. In particular, as /is varied, all the expected behaviors emerge: ultra-sensitive for /= 2, cooperative for /= 
1/2, anti-cooperative for / = — 1/2, non-cooperative for / = 0. 



This expression returns the average fraction of occupied sites corres- 
ponding to the equilibrium state for the system. We are now going to 
study separately the two cases of positive (/ > 0) and negative (/ < 0) 
cooperativity. 

The cooperative case: Chemical kinetics. When / > 0 interacting 
units tend to "imitate" each other. In this ferromagnetic context one 
can prove that the bipartite topology does not induce any qualitative 
effects: results are the same (under a proper rescaling) as for the Curie 
Weiss model; indeed, in this case one can think of bipartition as a 
particular dilution on the previous fully- connected scheme and we 
know that (pathological cases apart), dilution does not affect the 
physical scenario 35-38 . 

Differently from low-dimensional systems such as the linear Ising- 
chains, the Curie- Weiss model admits sharp (eventually discontinu- 
ous in the thermodynamic limit) transitions from an empty ((m A ) = 
(m B ) = 0) to a completely filled ((m A ) = (m B ) = 1) configuration as 
the field h is tuned. More precisely, eqs. 25 and 26 describe a trans- 
ition at a = 1 and such a transition is second order (Y changes 
continuously, but its derivative may diverge) when / is smaller than 
the critical value J c = 1, while it is first order (Y has a discontinuity) 
when / > J c . The latter case is remarkable as a discontinuous behavior 
is experimentally well evidenced and at the basis of the so-called 
ultra-sensitive chemical switches 39 . 

On the other hand, when / — » 0, the interaction term disappears 
and we expect to recover Michaelis-Menten kinetics. In fact, eq. 29 
can be rewritten as 



y(a;/) = 



l + ae 2/(2Y-ir 



(30) 



which, for / = 0, recovers the Michaelis-Menten equation Y (a) = a/ 
(1 + a) [Note: Note that we do not lose generality when obtaining eq. 
(7) and not Y (a) = a/(K + a.) (which is the more familiar MM 
expression) because we can rewrite the latter as Y (a) = iC^'a/Q + 
K _1 a) by shifting a — » a.K~ l and h — » [log(a/K)]/2.]. m this case there 
is no signature of phase transition as Y (a.) is continuous in any of its 
derivatives. 

In general, the Hill coefficient can be obtained as the slope of Y (a) 
in eq. 29 at the symmetric point Y = 1/2, namely 



«H = 



1 8Y 

y(i-y) aa 



r=i/2 : 



1 



(31) 



where the role of / is clear: a large /, i.e. / close to 1, implies a strong 
cooperativity and vice versa. 

One step forward, as the whole theory is now described through 
the functions appearing in the self-consistency, we can expand them 
obtaining polynomials at all the desired orders, more typical of the 
standard route of chemical kinetics. In particular, expanding eq. 30 at 
the first order in / we obtain 



(l-/)a + a 2 
l + 2(l-/)a + o! 2 ' 



(32) 



which is nothing but the Adair equation (eq. 2) as far as we set 
/= (l- VkM'kW 1 1) and we rescale a^a I \/kWkW. These 



results and, in particular, the expression in eq. 30 are shown in 
complete generality in figs. 1 and 2. 

The expression in eq. 30 can also be used to fit experimental data 
for saturation versus substrate concentration. Indeed, through an 
iterative fitting procedure, implied by the self-consistency nature of 
our theoretical expression, we can derive an estimate for the para- 
meter /and, from this, evaluate the Hill coefficient through eq. 31. As 
shown in fig. 3, fits are successful for several sets of experimental data, 
taken from different fields of biotechnology. The Hill coefficients 
derived in this way and the related estimates found in the literature 
are also in excellent agreement. 

Cooperative kinetics and cybernetics: Amplifiers and comparators. 
Having formalized cooperativity through statistical mechanics, we 
now want to perform a further translation in cybernetic terms. In 
particular, we focus on the electronic declination of cybernetics 
because this is probably the most practical and known branch. We 
separate the small coupling case (/ < J c , cooperative kinetics) from 
the strong coupling case (/ > J c , ultra-sensitive kinetics) and we 
mirror them to, respectively, the saturable operational amplifier 
and the analog-to-digital converter 40 . The plan is to compare the 
saturation curves (binding isotherms) in chemical kinetics with 
self-consistencies in statistical mechanics and transfer functions in 
electronics so to reach a unified description for these systems. 
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Figure 2 | Velocity of reactions versus substrate concentration or. These plots have been included to show the full agreement between our theoretical 
outcomes and the results presented in the celebrated paper by Levitzki and Koshland (see fig. 4 in Ref. 50). Different values of n H = 1/2, 1, 2, 
corresponding to /= — 1, 0, 1/2, are shown in different colors. Note that for this analysis there is complete proportionality between the reaction rates vand 
the saturation curves Y due to the law of mass action 50 . 



Before proceeding, we recall a few basic concepts. The core of 
electronics is the operational amplifier, namely a solid-state inte- 
grated circuit (transistor) which uses feed-back regulation to set its 
functions. In fig. 4 we show the easiest representation for operational 
amplifiers: there are two signal inputs (one positive received ( + ) and 
one negative received (— )), two voltage supplies (V sat , —V sat ) and an 
output (V out ). An ideal amplifier is the "linear" approximation of the 
saturable one and essentially assumes that the voltage at the input 
collectors (V sat and — V sat ) is always at the same value so that no 
current flows inside the transistor, namely, retaining the obvious 
symbols of fig. 4, i + = !_ = 0 40 . Obtaining its transfer function is 
straightforward as we can apply Kirchhoff law at the node 1 to delete 
the afferent currents, hence i 1 + i 2 + i- = 0. Then, assuming R l = 
1 Q (without loss of generality as only the ratio R 2 /Ri matters), in the 
previous equation we can pose i\ = — V_, i 2 = (V out — V-)/R 2 and i- 
= 0 (because the amplifier is ideal). We can further note that V_ = 
V + and V+ = V in so to rewrite Kirchhoff law as 0 = — V in + (V out — 
V in )/R 2 , by which the transfer function reads off as 

V mt = GV in = (l+R 2 )V in , (33) 

where G = 1 + R 2 is called "gain" of the amplifier. 

Therefore, as far as R 2 > 0, the gain is larger than one and the 
circuit is amplifying the input (R 2 < 0 is actually thermodynamically 
forbidden, suggesting that anti-cooperativity, from cybernetic per- 
spective, must be accounted by the inverter configuration 40 , see next 
section). 

Let us emphasize some structural analogies with ferromagnetic 
behaviors and cooperative kinetics. First, we notice that all these 
systems "saturate". Indeed, it is very intuitive to see that by applying 
a magnetic field h > 0 to a collection of spins, they will (at least 
partially, depending on the noise level) align with the field, resulting 
in (m(h)) > 0. However, once reached the critical value h such that 
m \ = 1 > an y further increase in the strength of the field (i.e. any 

h>h) will produce no variations in the output of the system as all the 
spins are already aligned. Similarly, in reaction kinetics, once all the 
ligands of a given protein have bound to the substrate, any further 
growth in the substrate concentration will produce no net effect on 
the system. In the same way, given an arbitrary operational amplifier 
supplied with V sat , then its output voltage will be a function of the 



input voltage V in . However, there exists a critical input Vj„ such that 
V ou t = Vsat> an d when input is larger than V,„ no further amplifica- 
tion is possible; the amplifier is then said to be "saturated". Of course, 
the sigmoidal shape of the hyperbolic tangent is not accounted by 
ideal amplifiers, yet for real amplifiers ± V sat are upper bounds for 
the growth, hence recovering the expected behavior, as shown in the 
plot of fig. 6 / [Note: One may notice that the outlined amplifier is 
linear, while in chemical kinetics usually slopes are sigmoidal on lin- 
log plots. This is only a technical point, and to obtain the logarithmic 
amplifier it is enough to substitute R 2 with a diode to use the expo- 
nential scale of the latter 40 .]. Moreover, we notice that the transfer 
function is an input/output relation, exactly as the equation for the 
order parameter m. In fact, the latter, for small values of the coupling 
/ (so to mirror ideal amplifier), can be written as (see eq. 17) 

(m)~(l+/)fe. (34) 

Thus, the external signal V in is replaced by the external field h, and 
the voltage V out is replaced by the magnetization (m). By comparing 
eq. 33 and eq. 34 we see that R 2 plays as /, and, consistently, if R 2 = 0 
the retroaction is lost (see fig. 3) and the gain is no longer possible. 
This is perfectly coherent with the statistical mechanics perspective, 
where, if / = 0, spins do not mutually interact and no feed-back is 
allowed. 

Analogously, in the chemical kinetics scenario, the Hill coefficient 
can be written as n H = 1/(1 — J) ~ 1 + /, in the limit of small / 
(namely for J < J c = 1, which is indeed the case under investigation). 
Therefore, again, we see that if / = 0 there is no amplification, and the 
kinetics returns the Michaelis-Menten scenario, while for positive / 
we obtain amplification and a cooperative behavior. This leads to the 
conceptual equivalence 

(^iV i+/ M g=i+/ ~t^/)' 

hence the Hill coefficient in chemical kinetics plays as the gain in 
electronics. This implicitly accounts for a quantitative comparison 
between amplification in electronics and in biological devices. 

One step forward, if / > J c the equation (m(h)) becomes discon- 
tinuous in statistical mechanics just like the corresponding (ultra- 
sensitive) saturation curve in chemical kinetics: The analogy with 
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data set J tih ± An// n'tf data set J n'tf 

VCH 0±0.6 1.0 ±0.7 0.97 CaMKII 1.1 ±0.3 4.4 ±0.6 
VC I-II 0.3 ±0.1 1.4 ±0.2 1.64 

CaMKII 0.67 ±0.03 3.0 ± 0.3 3 ±0.3 

Figure 3 | These plots show comparison between data from recent experiments (symbols) and best-fits through statistical mechanics (lines). Data refer 
to non-cooperative and positive-cooperative systems 22,51 (left panel) and an ultra-sensitive system 52 (right panel). For the latter we report two fits: Dashed 
line is the result obtained by constraining the system to be cooperative but not-ultra-sensitive (that is, / < 1 ) , while solid line is the best fit yielding / ~ 1.1, 
hence a "first order phase transition" in the language of statistical mechanics. The relative goodness of the fits are X 2 coop ~ 0.85 and X 2 ultm ~ 0.94, 
confirming an ultra-sensitive behavior. The tables in the bottom present the value of /derived from the best fit and the resulting n H ; the estimate of the Hill 
coefficient taken from the literature is also shown. 



cybernetics can still be pursued, but with analog- to-digital converters 
(ADCs), which are the corresponding limits of operational 
amplifiers. 

The ADC, roughly speaking a switch, takes a continuous, i.e. ana- 
logue, input and has discrete (dichotomic in its basic implementa- 
tion) states as outputs. The simplest ADC, namely flash converters, 
are built through cascades of voltage comparators 40 . A voltage com- 
parator is sketched in fig. 6h and it simply "compares" the incoming 
voltage values between the negative input and the positive one as 



follows: Let us use as the negative input the ground ( V = 0) as a 
reference value (to mirror one to one the equivalence with chemical 
kinetics or statistical mechanics we deal with only one input, namely 
the substrate concentration a. in the former and the magnetic field h 
in the latter). Then, if V in is positive the output will be V sat > 0, vice 
versa, if V in is negative, the output will be — V sat < 0 as reported in the 
plot in fig. 6m, representing the ADC transfer function. 

An ADC is simply an operational amplifier in an open loop (i.e. R 2 
= =°), hence its theoretical gain is infinite. Coherently, this corre- 




data set J n H ± An// n 



Azide 


3.0 


± 


0.1 


0.25 


± 


0.06 


0.23 


± 


0.06 


Thiocyanate 


2.5 


± 


0.3 


0.29 


± 


0.02 


0.28 


± 


0.06 


Cyanide 


3.0 


± 


0.1 


0.25 


± 


0.06 


0.22 


± 


0.05 


Imidazole 


2.7 


± 


0.4 


0.27 


± 


0.03 


0.25 


± 


0.05 



Figure 4 | Several sets of experimental data (symbols) 53 are fitted by eq. 30 with minus sign (solid line). The values of /corresponding to the best fits are 
shown in the table together with the related estimates for n H according to eq. 31. The estimates for n H obtained via standard Hill fit are also shown. 
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Figure 5 | This figure summarizes all the analogies described in the paper: In the first row, pictures of three biological systems exhibiting cooperativity, 
namely Mitogen-activated protein kinase 14 (positive cooperativity, panel (a), Ca 2+ calmodulin dependent protein kinases II (ultra-sensitive cooperativity, 
panel (b), and Synaptic Glutamate receptors (negative cooperativity, panel (c) are shown. The related saturation curves (binding isotherms) are shown in 
the second row (panels (d), (e) and (f), respectively), where symbols with the relative error-bars stand for real data taken from 41,52,54 respectively and lines 
are best fits performed through the analytical expression in eq. 28, obtained from statistical mechanics. The related best-fit parameters are / = 0.14, / = 
1.16, / = 0.29, respectively. Notice that in panel (d) it is possible to see clearly the "saturation" phenomenon as the first and the last experimental points 
are far from the linear fit (red line), while are perfectly accounted by the hyperbolic tangent predicted by statistical mechanics (green line), whose 
correspondence with saturation in electronics is represented in panel (1). Notice further that in panel (e), we compared the ultra-sensitive fit (solid line), 
with a simple cooperative fit (dashed line): at small substrate concentration the latter case does not match, within its variance, the data points (so 
accurately measured that error bars are not reported), while the former case is in perfect agreement with data points. In the third row we sketch the 
cybernetic counterparts, i.e., the operational amplifier (panel (g)), represented as an inverted flip-flop mirroring the symmetry by which we presented the 
statistical mechanics framework (the standard amplifier is shown in fig. 3), the analog-to-digital converter (panel (h)) and the flip-flop (panel (i)). The 
(theoretical) transfer functions corresponding to the circuits are finally shown in the fourth row (panels (1), (m) and (n), respectively) for visual 
comparison with the second one. 



sponds to values of / — » 1 that imply a theoretical divergence in the 
Hill coefficient, while, practically, reactions are referred to as "ultra- 
sensitive" already at n^^l. Consistently, as /— » 1 the curve (m(h)) 
starts to develop a discontinuity at h = 0 (see Fig. 6, panels d, e,fi, 
marking the onset of a first order phase transition. As a last remark, 
despite we are not analyzing these systems in the frequency domain 
in this first paper, we highlight that, when using time-dependent 
fields, for instance oscillatory input signals, full structural consist- 
ency is preserved as all these systems display hysteresis effects at high 
enough frequencies of the input signal. 

The anti-cooperative case: Chemical kinetics and cybernetics. We 

can now extend the previous scheme for the description of a 
negative-cooperative system, by simply taking a negative coupling 
/ < 0. Hence, eqs. 25 and 26 still hold and we can analogously 




out 



nctfe f 

Figure 6 | Schematic representation of an operational amplifier. 
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reconstruct Y (a) = ((m A ) + (m B ))/N versus a, whose theoretical 
outcomes are still shown in figs. 1 and 2, and fit them against 
experimental results as shown in the plots of fig. 5. 

Again, it is easy to check that there are two possible behaviors 
depending on the interaction strength /. If / < / c , the two partial 
fractions n A , ns are always equal, but when the interaction is larger 
than J a the two partial fractions are different in a region where the 
chemical potential log a. is around zero, as shown in fig. 1. In this 
region, due to the strong interaction and small chemical potential, it is 
more convenient for the system to fill sites on a subsystem and keep 
less molecules of ligands on the other subsystem. The critical value of 
the chemical potential | log a\ depends on the interaction strength: it 
vanishes when the average interaction equals ] a and grows from this 
value on. The region where the two fractions are different corresponds, 
in the magnetic models, to the anti- ferromagnetic phase. When / < /„ 
the binding isotherm, plotted as a function of the logarithm of con- 
centration, has a form resembling the Michaelis-Menten curve, even if 
anti-cooperativity is at work. Conversely, when / > J c , in the region 
around a. = 1 the curve has a concavity with an opposite sign with 
respect to the Michaelis-Menten one. In particular, there is a plateau 
around a = 1, which can be interpreted as the inhibition of the system, 
once it is half filled, towards further occupation. 

Finally, in order to complete our analogy to electronics, let us 
consider the simplest bistable flip-flop 40 , built through two saturable 
operational amplifiers as sketched in fig. 6 i, such that the output of 
one of the two amplifiers is used as the inverted input of the other 
amplifiers, tuned by a resistor. This configuration, encoded in stat- 
istical mechanics by negative couplings among groups, makes the 
amplifiers reciprocally inhibiting because (and indeed they are called 
"inverters" in this configuration), for instance, a large output from 
the first amplifier (say A) induces a fall in the second amplifier (say B) 
and vice versa. Since each amplifier pushes the other in the opposite 
state, there exist two stationary stable configurations (one amplifier 
with positive output and the other with negative output and vice 
versa). Thus, it is possible to assign a logical 0 (or 1) to one state 
and the other logical 1 (or 0) to the other state which can be regarded 
respectively as low concentration versus high concentration of bind 
ligands in chemical kinetics; negative or positive magnetization in 
ferromagnetic systems, low versus high output voltage of flip-flops in 
electronics. In this way, as the flip-flop can serve as an information 
storage device (in fact, the information (1 bit) is encoded by the 
output itself), the same feature holds also for the other systems. 
The behavior of the two flip-flop transfer functions (one for each 
inverter) are also shown in fig. 6 n, where the two (opposite) sigmoi- 
dal shapes are displayed versus the input voltage. 

Still in fig. 6, those behaviors are compared with experimental data 
from biochemical anti-cooperativity and their statistical-mechanics 
best fits with an overall remarkable agreement. 

Possible extensions: Heterogeneity and multiple binding sites. Another 
point worth of being highlighted is the number of potential and 
straightforward extensions included in the statistical-mechanics 
modelization. In fact, as the literature of mean-field statistical- 
mechanics model is huge, once self- consistencies are properly 
mapped into the saturation curves, one can perturb, generalize, or 
adjust the initial energy (cost) function and check the resulting 
effects. 

As an example, we discuss chemical heterogeneity, which has been 
shown by recent experiments 41,42 to play a crucial role in equilibrium 
reaction rates. To include this feature in our theory we can replace 
h ^2 m e 1- 1^ w i m yj j/lfCj with h t drawn from, e.g., a Gaussian 
distribution 



such that for a —* 1 homogeneous chemical kinetics is recovered, 
while for o->0we get standard Gaussian distribution Af[h, 1] for 
heterogeneity. 

We fitted data from 41 through the self-consistencies obtained by 
either fixing a = 1, or by taking a as a free parameter; the results 
obtained are in strong agreement with the original one. In particular, 
the authors in 41 found a ratio R between the "real" Hill coefficient 
(assuming heterogeneity) and the standard (homogeneous) one as R 
~ 0.53, while, theoretically we found _R ~ 0.57 (and a ~ 0.3). 

Further, we notice that n H grows with ja — lj, namely the higher 
the degree of inhomogeneity within the system and the smaller n H , in 
agreement with several recent experimental findings 42-43 . 

As a last example of possible extension, we discuss quickly also the 
multiple binding site case, which can be simply encoded, at least 
within the cooperative case, by considering an interaction energy 
of the P-spin type as 



£(MU,h) = 



_ . Ni N 2 Np 
i—ii=\ ' i[ i 2 i p 



(35) 



P(hj) cc exp 



{hi-hf 



2(1 -a) 



which results in multiple discontinuities for the binding isotherms 
as for instance happens when considering surfactants onto a poly- 
mer gel 44,45 , where the affinity of the surfactants to the gel is coo- 
peratively altered by a conformational change of the polymer 
chains (and actually these systems show hysteresis with respect to 
the surfactant concentration, which is another typical feature of 
"ferromagnetism") . 

Discussion 

In this work, we describe collective behaviors in chemical kinetics 
through mean-field statistical mechanics. Stimulated by the suc- 
cesses of the latter in formalizing classical cybernetic subjects, as 
neural networks in artificial intelligence 15,16,18 or NP-completeness 
problems in logic 46 " 48 , we successfully tested the statistical mechanics 
framework as a common language to read from a cybernetic per- 
spective chemical kinetic reactions, whose complex features are at the 
very basis of several biological devices. 

In particular, we introduced an elementary class of models able to 
mimic possibly heterogeneous systems covering all the main chem- 
ical kinetics behaviors, namely ultra-sensitive, cooperative, anti- 
cooperative and non-cooperative reactions. Predictions yielded by 
such theoretical frame have been tested for comparison with experi- 
mental data taken from biological systems (e.g. nervous system, 
plasma, bacteria), finding overall excellent agreement. Further- 
more, we showed that our analytical results recover all the standard 
chemical kinetics, e.g. Michaelis-Menten, Hills and Adair equations, 
as particular cases of this broader theory and confer to these a strong 
and simple physical background. Due to the presence of first order 
phase transition in statistical mechanics we offer a simple prescrip- 
tion to define a reaction as ultra-sensitive: Its best fit is achieved 
through a discontinuous function, whose extremization through 
other routes is not simple as e.g. least-squares can not be applied 
due to the discontinuity itself. 

It is worth noticing that, despite we developed mean field tech- 
niques, hence we neglected any spatial structure, we get a direct 
mapping between statistical-mechanics and chemical kinetics for- 
mulas, in such a way that we can derive from the former a simple 
estimate for the Hill coefficient, namely for "effective number" of 
interacting binding sites, in full agreement with experimental data 
and standard approaches. 

One step forward, toward a unifying cybernetic perspective, we 
described a conceptual and practical mapping between kinetics of 
ultra-sensitive, cooperative and anti-cooperative reactions, with the 
behavior of analog-to-digital converters, saturable amplifiers and 
flip-flops respectively, highlighting how statistical mechanics can 
act as a common language between electronics and biochemistry. 
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Remarkably, saturation curves in chemical kinetics mirror transfer 
functions of these three fundamental electronic devices which are the 
very bricks of robotics. 

The bridge built here inspires and makes feasible several chal- 
lenges and improvements in biotechnology research. For instance, 
we can now decompose complex reactions into a sequence of ele- 
mentary ones (modularity property 5 ) and map the latter into an 
ensemble of interacting spin systems to investigate potentially hid- 
den properties of the latter such as self-organization and computa- 
tional capabilities (as already done adopting spin-glass models of 
neural networks 49 ). Moreover, we can reach further insights in the 
development of better performing biological processing hardware, 
which are currently poorer than silico-made references. Indeed, from 
our equivalence between Hill and gain coefficients the more power of 
electronic devices is clear as G can range over several orders of 
magnitude, while in chemical kinetics Hill coefficients higher than 
n H ~ 10 are difficult to find. This, in turn, may contribute in devel- 
oping a biological amplification theory whose fruition is at the very 
basis of biological computations 6 " 8 . 

We believe that this is an important, intermediary, brick in the 
multidisciplinary research scaffold of biological complexity. 

Lastly, we remark that this is only a first step: Analyzing, within 
this perspective, more structured biological networks as for instance 
the cytokine one at extracellular level or the metabolic one at 
intracellular level is still an open point and requires extending the 
mean-field statistical mechanics of glassy systems (i.e. frustrated 
combinations of ferro and antiferro magnets), on which we plan to 
report soon. 
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